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Abstract. We study the primordial bispectrum of curvature perturbation in the uniform- 
density slicing generated by the interaction between the inflaton and isotropic background 
gauge fields. We derive the action up to cubic order in perturbation and take into account all 
the relevant effects in the leading order of slow-roll expansion. We first treat the quadratic 
vertices perturbatively and confirm the results of past studies, while identifying their regime 
of validity. We then extend the analysis to include the effect of the quadratic vertices at all 
orders by introducing exact linear mode functions, allowing us to make a reliable prediction 
long after horizon crossing where the features of both power spectrum and bispectrum are 
drastically different. It is shown that the spectra become constant and scale invariant in the 
limit of large e-folding, which implies the model can be consistent with the observational 
constraints regardless of the magnitude of the background gauge fields. It is found that 
depending on the period of inflation that falls into the observable window, the value of Jnl 
in the squeezed limit may well be within the reach of Planck. 
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1 Introduction 

The prediction on the primordial density fluctuation from inflation offers an exciting oppor- 
tunity to test the physics at high energy that is inaccessible for ground-based experiments. 
The advent of the Planck satellite, which is capable of improving the constraint on the three- 
point and higher correlations of the density perturbation at recombination by a factor of 10 
to 100, has prompted detailed theoretical investigations into the interaction of the inflatonfl]. 
So far, the efforts have been focused on the scalar self-interactions and interactions among 
multiple scalar fields. It is found that single scalar models with canonical kinetic term gener- 
ically predict an undetectable level of non-Gaussian signals [2] while a scalar with the DBI 
action or multi-scalar dynamics, such as hybrid inflation or curvaton scenarios, can lead to 
significant higher-order correlations [3-5]. These models being based on the string theories in 
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their origin of the inflaton, a detection of significant bispectrum or trispectrum can give us 
a clue for understanding the high-energy physics. 

In the context of unified theories of fundamental interactions, however, scalar fields can- 
not be the only ingredients of the universe. Most of the proposed theories such as superstring 
theories and M-theory rely on gauge symmetries, and gauge fields are indispensable to medi- 
ate interactions among the fields and in some cases to preserve supersymmetries. Even when 
they are absent in the fundamental Lagrangians, it is a generic prediction of dimensional 
reduction that a typical scalar field is coupled to some gauge fields [6, 7]. Earlier attempts 
to drive inflation with vector fields [8-11] turned out to be largely unsuccessful since one 
needs to abandon gauge symmetries, which results in the introduction of additional degrees 
of freedom and various instabilities [12-16]. More recently, interactions between the infla- 
ton and gauge fields, motivated by those unified theories of interactions, have been taken 
into account in the context of preheating[17-21]. In addition to interesting phenomenolo- 
gies including non-Gaussianity, primordial magnetic fields and gravitational waves[22-29], it 
was realized that the back reaction of gauge fields on the inflaton can effectively act as an 
extra friction term so that they slow down the rolling of the scalar field and help causing 
an accelerated expansion [30-36]. In fact, this back reaction can be so strong that it may 
generate a significant vacuum expectation value of the gauge fields and violate the isotropy 
of the universe. On the other hand, there has been a growing interest to maintain a small, 
but non-vanishing amplitude of classical gauge fields during inflation in order to explain 
the reported statistical anisotropy of Cosmic Microwave Background Radiation (CMBR) in 
WMAP 7- year [37-39]. By taking into account the aforementioned back reaction classically, 
the same types of scalar-gauge interactions arising from the high-energy particle theories have 
been found to enable acceleration of the cosmic expansion without requiring a sufficiently 
flat potential for the inflaton[40, 41]. This scenario turned out to be free of any classical 
instabilities or fine-tuning [42-4 7]. There have also been extensive studies on its potential im- 
prints on CMBR and it was revealed that even a very small amplitude of background energy 
density of the gauge fields could result in a significant statistical anisotropy in the curvature 
fluctuation [48-52]. While it implies such an anisotropic vacuum expectation value of gauge 
fields must be severely constrained, a recent study suggests that their effect on primordial 
bispectum is as drastic as its linear counterpart and the resulting non-Gaussianity may still 
be observable. In another recent development, it has been shown that multiple vector de- 
grees of freedom generically suppresses the residual anisotropy of the background space-time 
through a dynamical attr actor mechanism. In particular, when three or more gauge fields are 
coupled to a scalar field via a common gauge-kinetic function, the final state of the universe 
is completely isotropic regardless of initial conditions [53]. Such a circumstance may natu- 
rally be realised by non-Abelian gauge fields since the equal coupling is guaranteed by the 
symmetry [54]. There are other instances of isotropic inflation involving non-Abelian gauge 
fields which also exhibit similar attractor behaviours [55-60]. The linear perturbation of this 
isotropic inflation with background gauge fields has been studied, which has found that the 
primordial power spectrum is not strongly constrained by the current observations since the 
spectrum is perfectly isotropic and almost scale-invariant [61]. 

In the present paper, we investigate the second order perturbation of an isotropic uni- 
verse containing three U(l) gauge fields and a scalar inflaton. We compute the bispectrum 
of curvature perturbation by deploying the in-in formalism and compare the results with the 
corresponding work in an anisotropic background [52], which is expected to be qualitatively 
similar. There are theoretical, phenomenological, and technical reasons for this particular 
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model to be studied: 

1. The isotropy maintained by a triad configuration of gauge fields appears to be a generic 
feature of multiple vector degrees of freedom according to ref.[53]. This model serves 
as a prototype of the more complicated instances with non-Abelian gauge fields, for 
which similar features are expected. 

2. Because of its isotropy, the model cannot be effectively constrained by the power spec- 
trum. As we anticipate a strong signal in the bispectrum given the analogy to the 
anisotropic models, it is important to quantify it from a phenomenological point of 
view. 

3. While perturbation around anisotropic backgrounds is extremely involved, one can 
make a transparent perturbative expansion in the present isotropic model and identify 
all the relevant contributions. 

Besides, it is worth emphasising that the interactions under discussion frequently appear in 
supergravity theories that are low-energy effective theories of superstring theories and M- 
theory. It is therefore of great interest to study observational consequences of these models 
as they are beyond the reach of any ground-based experiments. 

Our results reproduce the previous studies when the e-folding number is relatively small 
while extending the analysis so that it is applicable to the period long after the horizon 
exit. A full perturbative expansion of the Lagrangian up to cubic order is carried out and 
several interaction terms that are not suppressed by any of the slow-roll parameters are 
identified. We first treat all the interaction terms, both quadratic and cubic, as perturbation 
and compute the three-point function for the curvature £. The amplitude is solely controlled 
by the parameter X 2 that represents the ratio of background energy density of the gauge 
fields to the scalar kinetic energy density. We explicitly show that the leading contribution 
comes from the vertex involving a scalar field and two gauge fields, which confirms the claim 
of ref.[52]. The three-point function scales as oc X 2 iV| where is the e-folding number 
after horizon exit for a mode with wavenumber k. The shape is local as has been shown 
in the previous studies. This results in a large /jvx when the modification to the power 
spectrum is assumed to be small. However, we find this conclusion is valid only if X 2 <C N^ 2 , 
which is not satisfactory for this isotropic model since X 2 is not necessarily that small in 
contrast to the anisotropic cases where this is required to keep the background anisotropy 
within the range allowed by the observations. A reason for the limited applicability is the 
quadratic vertices that generate an infinite number of Feynman diagrams in the perturbative 
expansion even at tree-level. In the second half of the paper, we take into account this fact by 
introducing exact linear mode functions. It turns out that one can solve the linear evolution 
equations analytically at superhorizon scales. By exploiting their general features, we shall 
prove that both power spectrum and bispectrum are convergent in the limit — > oo. Their 
final amplitudes are approximately determined by matching the superhorizon mode functions 
with the de-Sitter ones at horizon crossing. In order to obtain more quantitative estimates 
and handle the intermediate regime, we also solve the linear equations numerically from 
deep inside the horizon and use them in the integrand of the three-point correlators. We 
confirm the initial logarithmic behaviours in both power spectrum and bispectrum and their 
convergence at late times. It turns out that the time evolution f^i (squeezed) displays some 
interesting features. It first peaks at ~ 0.3X -1 where the peak value scales as X -1 , which 
means it may well be observable by Planck for a small X. Then /jvl monotonically decreases 
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and converges to a negative X-independent constant which appears to be —5/3. Therefore, 
this model is observationally viable regardless of the value of X as long as our observable 
length scales experienced a sufficient e-folding number after horizon exit. 

The paper is organized as follows. In the next section, after sketching the dynamics 
of the background evolution and introducing relevant parameters, we derive the perturbed 
Lagrangian up to cubic order. Section 3 gives the detailed procedure of computing the 
three-point function by perturbative expansion with respect to free de-Sitter mode functions. 
We calculate all the relevant contributions and identify the leading order term. Section 4 
discusses the importance of the deviation from the de-Sitter mode functions on superhorizon 
scales. By matching the solutions deep inside and far outside the horizon, we estimate the 
final amplitudes of the mode functions and the e-folding time when the analysis of section 
3 breaks down. In section 5, we confirm the analytical results and make the prediction on 
non-Gaussianity quantitative. Concluding remarks are given in section 6. 



2 Perturbative expansion up to cubic order 

Our model contains a scalar field and several gauge fields minimally coupled to gravity: 




( -^R - \d^d^ - V(<p) - iMlp^F^ 



R is the Ricci scalar curvature and F? v = d^A® — d u A^]a = 1,2,3 are three copies of 
U(l) gauge field strengths. These types of actions have been well studied in the context of 
magneto-genesis, preheating and anisotropic inflation. It has been realised that the coupling 
between the scalar field, which is identified to be the inflaton, and gauge fields enables an 
accelerated phase of expansion even with a relatively steep potential, as we will see later. We 
note that the energy density of the gauge fields stays constant in the first approximation in 
the inflating universe, violating the cosmic-no-hair conjecture. It has also been shown that 
the isotropic configuration of the gauge fields is a dynamical attractor of the system. Based 
on this result, we study the perturbation of this theory around the isotropic background with 
a non- vanishing triad of the gauge fields. 

We set 8nG = 1 and follow the ADM formalism [62] and parametrise the metric as 



t2 + N k N k 
Ni gij 



I -N z + N k N* Nj 

9^u ~ 



where 

N^g^Nj, g ik g kj 6' r 
The normalised extrinsic curvature of the constant time slice is given by 

E v = ~2 (sij - 2N m) 

and its intrinsic scalar curvature is 

3 R = {g ijM + feni^n,) (V V - g ij g kl ) ■ 

Electric fields are defined to be 

Ef = F a 0i . 
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2.1 Gravity and the scalar field 

The Einstein-Hilbert action in ADM formalism is given by 



2N v 3 



2^ , NJg 3 R 



E A \ + 



We assume the background is a flat Friedmann-Lemaitre-Robertson- Walker space-time 

ds 2 = a(rf) 2 {-drj 2 + 5 i jdx i dx j ) 
and write the perturbed metric components as 



N = a(l + 



a 2 /3i, gij = a 2 (5ij + 2^). 



When the problem concerns perturbation beyond linear order, one has to be careful in choos- 
ing the small quantities with respect to which the order of perturbation is determined. In the 
present case, we will solve the constraint equations so that <fi and are expressed in terms 
of jij and the other matter variables. Thus, their order in perturbative expansion is subject 
to the equations to be solved and we should distinguish different orders as 



Pi 



0(1) + 2^(2) + 



On the other hand, we avoid a similar expansion for 7^ since it will hardly appear in the 
following analysis as we are primarily working in the flat gauge, where the perturbation is 
set to be zero at each order by the choice of gauge. The only exception is the curvature 
perturbation on the uniform-density slice that will be expanded in terms of the dynamical 
variables in the flat gauge . As usual, the scalar-vector-tensor decomposition is made in order 
to decouple the linear order equations. It is defined by 



1 



-ifrSij + E tij + + -h 



S 



(n) 



0, 



ha = 0, 



h 



PI 
■ 0. 



B 



s 



(n) 



In the uniform-density gauge, we denote the curvature perturbation — £ = tp and expand it 
as 

C = C(i) + 2-C(2) + ■ ■ ■ ■ 

The 1 + 3 decomposition of the action for the scalar field is given by 



2N 



if' 2 - 2<f/<p ti N* + {<p,iN*) - N^g 



9 i:i V>,i<P,j + V((p) 



where primes denote derivatives with respect to the conformal time r]. We split cp into the 
background and the perturbation: 

if = If + 7T. 

7r will be treated as the dynamical variable in terms of which the perturbative expansion is 
defined so that we do not need to expand it further. 
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2.2 Gauge field perturbations 

The Maxwell Lagrangian in ADM formalism reads 

Cm = ^f 2 9 13 (E% + F a lk N k ) (E a 3 + F^N 1 ) - ^-fY k g jl F%F%. 
The perturbative expansion of the vector potentials yields 

A a = a a , At = A(r,)5 a l + X 



a 



where the background quantity A(rj) behaves effectively as a second scalar field. The back- 
ground values of Aq are taken to be zero by a gauge choice. As in the gravity sector, we 
should in principle distinguish the different orders of perturbation for the variables that are 
expanded in terms of the dynamical ones. However, after the adoption of flat slicing and 
£7(1) gauge fixing, we are left only with the dynamical variables from this sector. Hence we 
suppress this distinction and the scalar- vector-tensor decomposition is carried out as follows: 

c a = M,a + "a, X a i = a5 ai + oi + e aij (tj + Xj) + + u ai , 

= Xij = Ki t i = 0, Ula = 0, = 0- 

2.3 Background dynamics and parameters 

Before going into the perturbative analysis, we briefly review the background evolution of 
the system and identify the relevant parameters. The Maxwell's equation can be trivially 
integrated to give 

A'- — 

where c is an integration constant. As usual, we introduce the "slow-roll" parameters 

« = » = & (* = 7)- < 21 > 

which characterize the evolution of the scale factor a(rj). The Raychaudhuri equation 

2H' + H 2 = -hp' 2 ~ + a2 V ( 2 - 2 ) 

tells that the potential energy has to dominate over the scalar kinetic energy and the energy 
of gauge fields in order to have an accelerated expansion. This suggests the introduction of 
another parameter 

(2.3) 



v 2H 2 

which controls the evolution of the inflaton. Combined with the Friedmann equation 
one derives 

c 2 
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which is the representative of the energy density for the gauge fields. Note that e v < en 
where equality holds when the gauge fields vanish. Since this deviation from the single-scalar 
inflation plays a central role, we define the parameter 



Z=.p— Sp, (2.5) 

V e <p 

which measures the ratio between the energy density of the gauge fields and the scalar kinetic 
energy. We note that X does not have to be small as far as the background dynamics and 
the power spectrum are concerned. Without loss of generality, we can assume <p' > and 
use (f> = ■y / 2e tp / H. Now the equation of motion for (p gives 

v" = ^{2-2e H + i^)U 2 (2.6) 

where 

In principle, this quantity does not have to be small as long as ijh <C 1, but we do assume it 
is in order to control the perturbative expansion. Now by differentiating 



c 



2 



P 

one obtains 



J ~ 



P 2(e H 
and using the equation of motion for the scalar field yields 



a 2 Vu> If 2 3 1 \ . . 

6e H - 3e H + e H e<p + ^e H r}H - J ■ I 2 - 9 ) 



The first expression tells that the slope of f{ip) must be steep in order to maintain the 
amplitude of gauge fields during inflation. The second implies the gradient of potential is 
not necessarily small if <C e#, or equivalently, if X S> 1. The reason is that the slow-roll 
of the inflaton can be achieved by transferring the scalar kinetic energy to the gauge fields 
through the coupling /(</?). It later turns out that the perturbative approach breaks down 
when X > 1 anyway, so we assume that X < 1, where the usual intuition from single-scalar 
model works well. The higher-order derivatives of V and / take complicated forms in general, 
but assuming the constancy of i]H,ip and keeping only the leading order terms in the small 
parameters, we obtain 

-(4e H -2 VH + Vip ), (2.10) 



U 2 2e 



v 



fl2 ^w 3 €H ( 



V? 

and 



3 

(/ ), w 8 (/ ),<p tp< p 



p ^ p w ■ (2 ' 12) 

Finally, we emphasize that this regime of accelerated expansion aided by gauge fields is a 
dynamical attractor for a wide range of potential and coupling. The readers are referred to 
ref.[53]. 
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2.4 The cubic action for scalar perturbations 

Since the curvature perturbation does not receive any contribution from vector or tensor 
modes at linear order, we can eliminate them from the tree-level calculations of three-point 
correlation functions arising from cubic interactions. Since their power spectra are known 
to be small, the higher order contribution is expected to be negligible. Hence, we focus on 
scalar perturbations hereafter. 

For the scalar perturbation, the gauge field variables are given by 

X a i = Ot5 ai + ai + e ai jTj. 

We use the U(l) gauge freedom to set 6' + // = 0. It follows that 

^^ij = &,j3ia Qjifija ^aikT^kj ^ajkT,ki^ 

5Ef = a' + e aij T'j, 

where 5 indicates the perturbation of the following variables. For the metric, we adopt the 
flat slicing tfj = E = 0. Focusing on scalar modes, we can completely ignore ^j. The 



gravitational action drastically simplifies up to cubic order to become 

a~ 2 C g = (l-0(i)) 
The scalar part is the same as the standard: 

£y> = (1 - 0(1)) Q^^l) - 4>{l) (<£V + aV^Tr) - (p'TT,iBf + \' 2 - ~ ^^W^j 

-a2V tp 4> 2 1 -)ir - tt'tt^B^P - ^a 2 V wip ir 3 . 
After some straightforward algebra, one obtains the gauge Lagrangian as 

C M = (1 - (1) )4% - 2fa,Bf (a> + ^^vr) + f 2 B f (e ljk a^ k - r y rj - t't n ) 



+ c 2 (/ 2 W ^3 + HQw^, + { f)^ (la' 2 - a % a 



where 

f2\ \ 0^2/^2 



A2) 3c 2 2 ( , c(/% \ (1) 3c 2 (/ 2 ), w 2 3c(/ 2 ),,, , 

+— (3a' 2 - 2a ;a; + 2r / i r / i - t^-t^ - r^r . 
Therefore, the total Lagrangian up to cubic order is written as 

C = (1 - (1) )£( 2 ) - aV,4)i - aVTT^W _ 2 / 2 (a' + ^p^) a,B 

C 2 {f 2 ) ,ipipip 1 4 T/ \ 3 3c(/ 2 ) jW q / | //?2\ „ /3 /2 



/ 2 



(1) 



+ ( J U 4/ r - + ^^a' + (A^tt ^ - a >4 a,, j (2.13) 
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where the quadratic Lagrangian is given by 

.,2 



-<ft*,#f - + V - l w - \a 2 V^ 2 (2.14) 

, 3c 2 (/ 2 ) w 2 , 3c(/ 2 ) ^ , f 2 / /2 \ 

+ 4a 2 /4 ^ + a 2 /2 a7r + ^2( 3a -2a,iO,<) 

It should be mentioned that we dropped the terms involving <^ 2 ) an d -B^ 2 ) from the beginning 
since they multiply the background and linear-order constraint equations, which would be 
automatically satisfied in our formulation. 

2.5 Solving the linear constraints 

Using the background equations and parameters, the quadratic Lagrangian can be rewritten 



MS 



a 



2 £ (2) = -a 2 V4> 2 {l) - U<P {1) f 2V 2 B (1) + T^vr' + ^-^la'^j + ^H 2 (1) vr 
+HB^ f + ^V 2 «j + V 2 1 • 



— 7T .;7T j 
2 ,» 



( 3a ' 2 ~ 2a ,i«,i + 2r ,i T ,i - T ,ij T ,ij ~ T ,U T ,jj) 



where we discarded the surface term and defined 

1 



66^,(1 + 2X 2 ) - 6e 2 H + 4e^e v + 36/^ - 2e ip rj lp ) . 
Varying B^ determines 0m as 

= V2 + ^f Xa ) • (2 - 15) 



Using this, variation of 0m leads to 



^-V 2 B^ = -V + ( 6 J 2 + ^X 2 (5 - 6X 2 ) + \ m (l - X 2 ) - rj^j U>* (2. 



16) 



— ■=-! (3a' + (6 — 3en — e ip )'Hoi\ . 
V2a 

These relations will be substituted into the Lagrangian derived in the previous subsection 
and the curvature perturbation in the uniform-density gauge introduced in the next. 
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2.6 Curvature of the uniform-density surface 

For the purpose of quantum field theory calculations in the multi field dynamics, the most 
convenient gauge is the fiat gauge where 7^ = hij [3]. However, the observationally relevant 
quantity is the curvature perturbation in comoving gauge 1Z C that coincides with the curvature 
in uniform density gauge £ beyond the horizon scale. The latter is more often picked up 
as is done here since it possesses a desirable mathematical property. Hence, we need the 
transformation law between fiat gauge and uniform density gauge, which we cite from ref. [63] 
as 

-C(i)=n^ (2.17) 

and 

" C(2) = j \SP(2) ~ ^T^(l)) - + Jv- 3 Stf,y (2.18) 

where the right-hand-sides are evaluated in the flat gauge. We defined the perturbative 
expansion of the energy density 

p = p + 5p {1) + -6p(2) + ••• 

and a quadratic expression 

~p { S PW,i B f + S PW,j B ?) - p2 5 Pm 5 P0),i- 
The background sound speed c 2 . in the present setting is 

2 P' 2 1 

c s = — = -l + -e H + -m 
p' 3 3 

with p being the background pressure. At the linear order, the energy density in the flat-gauge 
is neatly written as 

a 2 Sp {1) = -6n 2 (f) {1) - 2nV 2 B^\ 
The background energy density satisfies 

thus 



p = -6e H 



a 2 " 



Therefore, the first order curvature perturbation is given by 

C(i) = -^(3^ (1) +V 2 i?( 1 )). (2.19) 
The second order part will be discussed later. 
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3 Analytical estimate of the bispectrum in the limit of small X 



In this section, we apply the standard methods of the in-in formalism to the Lagrangian 
obtained in the previous section. In order to render the problem tractable, we keep only the 
leading order contributions in the small parameters €H,(p,VH,tp- An interesting point is that 
even in the limit of de-Sitter space-time, the key parameter I does not necessarily vanish. 
Using the eqs. (2.1) - (2.12) and substituting (2.15) and (2.16), the cubic Lagrangian (2.13) 
becomes 

a~ 2 £ = I (vr' 2 - (Vvr) 2 ) + £ (3a' 2 - 2(Va) 2 ) + ^vr 2 + ^4ra' (3.1) 
2 zcr rj z rj a 

^rvr 3 -vrV - J — ^tt (3c/ 2 - 2 Vaf . 

e,on z ^fe~ n n a \ e, a a 2 



We dropped r for a reason that becomes clear soon. We assume the background is close to 
de-Sitter, which implies 

a = -^> f = fw 2 (3-2) 

and discarded all the terms higher order in slow roll. One can rescale a to set fo = 1 without 
loss of generality. We further demand I < 1 since we would like to treat all but kinetic terms 
perturbatively. The factors of y^ -1 appearing in the cubic terms might look worrying for 
validity of the perturbative approach. But when the action is written in terms of £, they are 
of the same order as the quadratic kinetic terms and the perturbative expansion should be 
marginally applicable. The following analysis is expected to be valid for en Cl< 1. We 
are concerned with the three-point correlation function of the curvature perturbation 

which is obtained from (2.19) with (2.15) and (2.16), neglecting all the higher order terms in 
slow-roll. The absence of r at this linear order justifies its omission from the Lagrangian. 

3.1 Notations 

We take the free massless part of the action to be the background and treat all the other 
terms perturbatively. In the interaction picture, we set 

Tr/fox) = J (03^|= (M^lWe^ + ul^ale-^) , (3.4) 

= £ / wf^m <S (7?)6ke *' x + ^ e ~ ikx ) • < 3 - 5 ) 

where the de-Sitter mode functions are defined to be 

Mv) =(i~ kr])e- ikri , 

v k {v) = {cskrj - i)e- iCakr > 
The interaction Hamiltonian is given by 

r2 
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where 



6V21 



Hj = ^ 

H f = - JpL J d 3 w Trja'j, Hf = 3 J d 3 w irja'f . 

The term with higher spatial derivatives has been omitted. We often drop the subscript I. 
Note that it was claimed in ref.[52] that Hf gives the leading contribution to the bispectrum. 
We shall explicitly confirm that it is the case as long as we remain within the regime of validity 
for the perturbative treatment of the quadratic vertices (i.e. Hj and the mass term for ir). 
We are going to compute the three point correlation function in Fourier space defined 

by 

<Cfa,x)Cfa,y)Cfo«)) =/// (2^)3 H ^ (OA^ 

We often abbreviate it as (C k )- 

For later use, the following mode functions are also introduced: 

e(r,x)= 7 r' + ^ = / ( ^)3 ^ {^We^ + gt( V )ate-^ 

g k (r,) = [i(3 + k 2 rf) - 3k V ] e~ ikr >, 
h k ( V ) = [i(3-c 2 s k 2 V 2 )-3c s k V ] e- ic ° kr >. 

It will be useful later to note that 

g* k (v)u k (v) = (3 + k 2 r](r] + 3fj) - 3ik(r) - fj) + ik 3 rj 2 fj) e ik ^-^ , (3.6) 
K(v)v k (v) = ("3 + c 2 s k 2 ri(r] - 3fj) + 3ic s k(i] - fj) + ic 3 s k 3 T] 2 fj) e 1 ^'^ . (3.7) 

3.2 The outline of the calculation 

The three point function for Q is written in terms of ^ and a' as 

H (C(r ? ,x)C(r ? ,y)C(r/,z)) = rf{^, x)£fa, yMv, z)> (3-8) 

3ifZ?7 6 



((C(^,x)^(77,y)a'(r?,z)) +2 perms) 



H ^ ((C(»7 ) x)a / (77,y)a / (7? ) z)) +2 perms) . 



Despite the appearance of lower powers of I in the Lagrangian, the leading order contribution 
to (( 3 ) turns out to be quadratic. 1 Then, the tt 2 term in the interaction Hamiltonian is clearly 



1 This is true for tree-level calculations, but loop contributions may contain terms without any factor of T. 
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irrelevant. However, we do have to keep Hj since it affects, for example (£ 3 ) with Hf at this 
order. More specifically, (£ 3 ) can be written as follows: 



(e) = ij d m ([Hf( m ),ej}} (3.9) 

/v rm 
d m J d m ([Hf(r, 2 ),[H](m),ei]]) (3.10) 

/T] fill 
dm J d m ([H«( V2 ),[H?(m),tf}}) (3.ii) 

/v rm rm 
dm J d m J d m {[Hf{ m \[Hj{m)A H K^UW ( 3 - 12 ) 

/v rm rm 
dm J dm J d m ([H]( V3 ),[H?( m ),[H<>(m),ei})]) (3.13) 

/m rm 
dm J d m ([Hj( m ), [HKm), [#f faiU?]]]) + o(i 3 ). (3.14) 

In a similar way, (£ 2 a') contains the quadratic term in the Hamiltonian given as 

<eV) = i J v dm([H?(m),£<x'i]) (3.15) 

rv rm 

-J dm J d V2 ([H?( m ),[H' 1 I (m),eiu' I ]]) (3.16) 

/v rm 
dm J dr ?2 ([^( ??2 ),[ J fff(r ?1 ),e / M]]) + 0(X 2 ). (3.17) 

Finally, (£c/ 2 ) receives no contribution from the quadratic interaction and becomes 

(Ca' 2 ) = i ^ dm([H<f(m),Zia?}) + 0(l). (3.18) 

Hence, we have to compute 10 distinct integrations to fully work out (£ 3 ) at the quadratic 
order in X 2 . We focus on the superhorizon limit of the spectrum, i.e. —kiT] <C for all 
i = 1,2,3. 

3.3 1-vertex contributions 

Let us start from the single integration (3.9). One can easily go to Fourier space and derive 

. 3 . [2 4Z 2 H 6 r dm rl * * * * * * \ 

= *V e^ ^V 8fc 3 fc 3 fc 3 J ~4 { u kiUk 2 u k3 g kl g k2 g k:j - g kl gk 2 gk 3 u kl u k2 u k3 ) 

2 6X 2 i"/^ /" r ' dt]i 

At a first glance, the integral looks divergent as rj — > even if the factor of if 1 in (3.8) is 
taken into account. However, it is not the case since we have 

3 « (v)g* k2 iv)g* k3 (vHi iv)u k2 {r])u ka (rf)) = 9 (kf + k 3 2 + fc|) rf + 0{rf). 
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Therefore, an integration by parts gives 

Ai= [ (i)4 (i?i)«fc 3 

= -3(fc? + kl + k 3 3 ) + 0( V 2 ) + \ f %Q (iii (u kl u k2 u k3 )') . 

The same type of cancelation of power holds for the remaining integrals. It is also helped by 
the fact that 

u' k (rj) = ik 2 V e- ikr >, 

which is a manifestation of the constancy after horizon exit of the de-Sitter mode functions. 
Repeating another integration by parts, we are left with 

Ai = 6(kf + kl + kl) + 0(r] 2 ) 
The third line is finite. The leading contribution is logarithmically divergent in r/ as 

.^ 3 . Him 4 ! 2 (kf + kl + kl) ndm 

(CD - ^- ^ k3 J J f cos [ {kl + k 2 + .3) (v - m)] • 

The cross correlations (3.15) and (3.18) are in principle similar. We have 

24X 1 

= / ^rh kl (r] 1 )u k2 (i]i)u ks (r]i). 

J m 

and 

/2 6 H 2 1 
-Z9^3 36c 6 fc 3 fe 3 S {9* kl (vMM h k(v) H ki(v)) , 
r • 1 S 2 3 

P dm 

HkAv) = / —u kl (rn)h k2 (m)h k3 (vi)- 

J m 

In carrying out these integrals, it is useful to note that 

?** = (?)'■ (3 ' 19) 

It will be later useful to derive the explicit forms of T k and H k . Straightforward integrations 
by parts lead to 



^ c -i(c s k 1 +k 2 +k 3 )ri 1 



T kl (v) = -{kl + kl) J' 
' i c s ki + k 2 + k 3 

j] 3 rj 2 r] ' c s ki + k 2 + h 



+ fj__ c sh + k 2 + k 3 _ . Cskx (k 2 + fc 3 ) - fc|_+ k 2 k 3 - fc| + c 2 s kfk 2 k 3 ^ 
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and 

H k (77) = -3k 3 r dVl c ~ i(Cakl+eek3+k3)ni + i c i k ik%k 3 -i{c s ki+c s k2+k 3 )r) 

3 J m Csh + c s k 2 + k 3 

| c 3 s hk 2 (c|fcifc2(fci + fc 2 ) + c,fc 3 (3fcf + 8fcifc 2 + 3fc 2 ) + 3(fci + fc 2 )fcg) c _ i(Cakl+Csk2+k3)TI 

(c s ki + c s & 2 + k 3 ) 2 

+ /3£ _ 3(c s fc 1 + c s fc 2 + fc 3 ) _ , 3 (c^ifc 2 + c,fc 3 (/ci + fcg) - fcf) \ e - i{Cs k 1+Cs k 2 +k 3 ) v _ 
\rj 3 r] 2 77 y 

All the terms with negative powers of 77 cancel when taking the imaginary parts, due to the 
rapid decay of the imaginary part of the propagators beyond the Hubble horizon; 

9 (gl(v)Mv)) = k 3 V 3 + 0{rf) and 9 (h* k ( V )v k ( V )) = c 3 s k 3 V 3 + 0{rf>). 

The end results are again logarithmic dependences on 77; 

r /v 27H 3 l(kf + k%) Pd m 

-> ~ ^=~6^3fcpp3 _/ — COS + fc 2 + C^ 3 ) (7/ - 7?l)] , 

r , v 27tf 2 fc? PdT/i 

(Cfci«fc 2 «fc,) -> „ S= o c, - TTsTs / — cos l(*l + c ^2 + c s /c 3 ) (77 - r/i)j . 



2 ^/2i^r7 9 cf fcf fc| fc| J 771 
3.4 2-vertex contributions 

Let us start from the term (3.10). By definition, the connected tree level contribution is 
given by 

-(H B ( m )Hi( m )e) = j d 3 w 2 ^- j d 3 Wl (^(w 2 ) 2 a'(w 2 )vr(w 1 )a'(w 1 )e(x)e(y)e(z)) 

[~2 144X 2 f f 

= J- jp~ J J d 3 w 2 d 3 w 1 (a'(w 2 )a'{w 1 )) ((vr(w 2 )^(y)(7r(w 2 )e(z))(vr(wi)^(x)) + 2 perms) . 

Given 

(2tt)3 2fcgiy "* 1 {VlUl ^) e ~ aCr(X ' Wl) » 
(a'(7 ?2 ,w 2 )a'(7 ?1 ,w 1 )) = I ^fl^ WW^^. 
we can write down its Fourier transform as 

-{H B {ri 2 )H q {rn)il) = J— 3 f5 3 ^ 3 ~ g* kl (r])g k2 (77)^ (77)^ (771) fcl ^ nfefaKaM ^ 

+2 pemrs. 
Taking the commutators, it becomes 



o 19 /T^T^ 

V e ¥> c s K l ^ 2 "-3 '/ 



xS ( fffc 2 ( r ?)fffc 3 fo) kl T^- u k2 {m)u k:i (772) ^ 4^ ) f 2 perms. 
^1 % 
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Now our task is to carry out the integral 

A 2a = f {ai^kAm)) 3 {9k 2 9*k 3 h ii(vi) T ki(m)) ■ 

J Vi 

Substituting the integrated expression for T^r/i), the single integral arising from its second 
line gives at most oc ln(— rj) in the limit rj — > since any power divergence disappears after 
taking the imaginary part as was demonstrated for the 1-vertex cases. 2 We then only have 
to check if the remaining term gives rise to similar logarithmic contributions. Evaluating the 
double integral, we obtain 



u kl + ikWe~ lklV 



9 ( 9 k2 9 ka v kl — e 



p -i(c s ki+k2+k 3 )yi 



>w- «.-•])»(«;,«— 

All the integrations are at most of order ln(— i]) as rj — )■ except for the first line whose 
leading term yields 

-!- / -i- cos(A;i (77 - 771)) cos (c a fci(?7i -772) + (k 2 + k 3 )(r)- 772)) 

and behaves as (hi(— rj)) 2 . The leading order behaviors for the other contributions are essen- 
tially the same. For (3.11), the amplitude reads 

-([H«(r, 2 ),[H B ( Vl Ul])) = ~\ - M Su3^ (V)9* k3 (??K 2 faKa fa)) 



rc l ^2 ^3 '/ 



xS UM^UkM^) +2 perms. 
V ^1 V 2 / 



Defining 



^fci(m) = / ^ru kl (r] 2 )h kl (r]2) 

J m 



2 In fact, one can directly show there is no power divergence by Taylor-expanding the propagators (3.6) and 
(3.7). Namely, the integrand does not contain any power of 77, 771 and r\i lower than —1 when the exponential 
is written as a power series. Since the calculation relies on Mathematica, however, we stick to the integrations 
by parts here. 
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we retain only the most divergent term to obtain 
J Vi 

~ -k\ f ^ (gi 2 gi z [kle- ik ^u k3 + kle~^u k2 ] ) 9 (ig* kl v* kl f" ^V^d+^A 
~ 27k\{kl + kl) f ^ r ^ cos ((A* + k 2 )(v - m)) cos (*i(i7 - m) + c s fci(m - %)) ■ 

J vi J m 

For (3.16), 

-( J f/ c ( ??2 )^( ??1 )^(x)e(y)a / (z)) = -^-^ jj d 3 w 1 d 3 w 2 (a'(w 2 )a'(w 1 ))(a'(w 2 )a'(z)) 

x ((vr(w 2 )C(y))^(w 1 )e(x))(7r(w 2 )C(x))(vr(w 1 )e(y))) , 

and in Fourier space, it becomes 

-(H C ( V2 )H«( Vl )Z kl Z k2 a' k3 ) = - ^« nW 9lM9lMKM 

Z V tl f 7 c s' i l'>'2 ft '3 



*"fci(*ft) 4 ^i feVfa fe 4 + 1 2 ■ 



Then, we derive 



n EJ3-7" 

-([i^faM^U^^]] = ^, 6j , 6 ,3. 3. ^(^(^ fa)) 

V fc ¥> t s' t l fi '2 /i 3 7 



^ 2 fa^(^)^^^i(^)^ 2 fe)^^) + (1 O 2). 



As before, integration goes as 



&2a = I (g kl Uki(m)) % (gl 2 hl 3 hl 1 (r) 1 )H k2 (i 1 i)) 

7i 

3kl f ^ (g* kl u kl ) 9 (g% 2 h% 3 h%\ f 1 ^ e -^k 1+ k 2+Cs ks) V2 
i dr) 



/v dm f Vl dr]2 
/ cos (ki(rj - rji)) cos (c s ki(iji - r/ 2 ) + (k 2 + c s k 3 )(rj - rj 2 )) 

Finally, (3.17) yields 

-(m m ),[H G ( m U kl tk 2 * k3 ]] = ' 6 f 6 3 ,l 3,6 (& (vM 3 (vW 2 (m ) 

(glMKM) 1 ^^^^ + (i o 2). 
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i{l+c s )kirj2 



The leading order contribution is 

®2b = [ ^t-Q {g*k 2 h t 3 u k 2 (m)hk 3 (m)) 3 (s^/t^M^kM) 

~ -k\kl r ^ (g* k2 h* k3 e- ik ^h k3 ( Vl )) 9 (ig* kl v* kl ( Vl ) H 

/ cos ((k 2 + c s k 3 )(r] - rn)) cos (cshirn - T] 2 ) + h(r] - T] 2 )) ■ 

3.5 3-vertex contributions 

We saw that the 1-vertex terms that involve only single time integrals resulted in oc ln(— 77) 
while the leading contributions from 2-vertex terms come from double integrals and propor- 
tional to (ln(— rf)) 2 . Hence, one expects that 3-vertex contributions behave like (ln(— r/)) 3 and 
dominate the tree-level amplitude at the order X 2 . This was also the result of the ref.[52]. 
We explicitly prove it and derive the coefficients in front. The principle of the calculations is 
the same as the previous sections although the algebra gets increasingly complicated. First 
of all, we write 

-i( J ff c (r ?3 )^(r ?2 )^( ??1 )^x)C(y)e(z)) = -i^Jl^- fff d 3 wfw 2 d 3 w 3 
x(a'(w 3 )a'(w2)}(a'(w 3 )a'(wi)) ((7r(w 3 )£(z))(7r(w 2 )£(y))(7r(wi)£(x)) + 5 perms) 

and 

3H 4 1 2 1 



(H c ( m )H\r, 2 )H^ m )^ k2 ^ 3 } 



^4 RT3TEU9k 1 {v)9k 2 {.'n)9k 3 {'n) 
26^3 (7717/2%) * cgfcJA§fc| 

xukAvi)KAvi) u k 2 (V2)hl 2 (m)hk 1 (ri3)hk2(V3)uk 3 (V3) + 5 perms. 



In the end, our integrand is 

{gl^kAm)) Q {gl 2 u k 2 {m)) 9 {gk 3 KS^) h t 2 {^) h kAm)K 2 {m)uk 3 {m)) + 5 perms. 

Anticipating the cancellation of terms with negative powers of 77, we seek the expected 
(ln(— 77) ) 3 contribution. It can only come from 



A 



3a 



(P*>fci(m)) / (5^^2(772)) 9 {9k 3 h *ki(vi)K 2 (V2)H ks (v2)) 
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where we integrated by parts for 772 • We perform another integration by parts with 771 as 
follows: 

^ 9 ( g * ( V 1M1\ ' v * {m) f V2 dsm e _ l{Cakl+Csk2+Csk3)m 



111 J — J 773 



d m %(g* kl u kl (r,i)) r ^(g* k2 e- 1 



Vl V ' V J V3 

Only the second term can give rise to the sought dependence on 77. We derive 

,3^3 r ^ (n* P -ikim\ p ^HIq (gi 2 e-^w 

-i(c s ki+c 3 k 2 +k 3 )r) 3 



xS ( 9ty kl (vi)vt 2 (m) J — e" 

/v d m rvi d v 

— cos (A* (77 -771)) / — cos (£2(77-772)) 

r n dm 

x / cos (c s fci (771 -773) + c s £2(772 -773) + £3(77-7/3)) . 

J V3 



For (3.13), we have 

-i(^(77 3 )if c (772)^(77 1 K(x)e(y)e(z)) = -avt^^r III d 3 Wl d 3 w 2 d 3 w 3 
x(a / (w 3 )a'(w2))(a / (w 2 )a / (wi)) ((7r(w 3 )f (z))(7r(w 2 )£(y))(7r(wi)£(x)) + 5 perms) 

and 

' 2 12H 4 1 2 1 



([H0( m ),[H c ( m ),[Hi( m ),a]]] 



^ v 3 (mmmf c6 s k i k 2 k l 

{gt^kAvi)) $ (Sfca^fciM^M^ite)) 9 {gi 3 K 3 (m)uk 3 {m)hk 3 {m)) + 5 perms. 
The integration results in 

^3fo = / — ^(fffc^fciCm)) / —^{9k 2 K 1 (Vi)uk 2 (m)hkAV2))^{gl 3 hl 3 (i 1 2)F k3 {r]2)) 



kl f ^3 (^ e - lfcl?n ) (^ 2 ^ 1 (m)e-*' ?2 ^ 1 (77 2 ) 



- £?£| 

^ d m 



/>7 c / 7 y 1 rvi 
cos (£1(77 - 771)) / cos (c s £1(7/1 - 772) + £2(7/ - m) ) 

— cos (£3(77 - 7/3) + c s £ 3 (t72 - 773)) . 
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For (3.14), we have 

-i(H^ m )H^ri 2 )H c ( Vl )^ny)^)) = -* 3 % 4 y~f- /// d s Wl d 3 w 2 d 3 w 3 
x(a'(w3)a'(wi))(a'(w 2 )a'(wi)) ((7r(w 3 )£(z))(7r(w 2 )£(y))(7r(wi)£(x)) + 5 perms) , 

and 

-i([H«( m ), [ H \ m ), [H c ( m ),e k ]}} = _ /1 3 12g4x2 1 

{g* kl u kl (m)) ^ {9k2 h tx('ni)uk2{m)hk 1 {m)) ^ {gk 3 h *k 3 (mhk 3 (m)hk 3 (v3)) + 5 perms. 

Similar to the other two, we obtain 



A 2c = / {g^Am)) — 9 (^^M^M^M) Q {g k3 h* k3 (vi) F k 3 {v2)) 

J Vi J V 2 

~ -klkl [' ^ (g* kl u kl (Vi)) ^ (9*k 2 h* kl ( Vl )e- ik ^v kl ( V2 )) 



?7i 

-i(l+c s )fc3r;3 



/ v drji f Vl dr}2 
cos(fei(7/ - 7?i)) / cos (c s ki(rji - rj 2 ) + k 2 (r] - m)) 

— cos (k 3 (r] - 773) + c s fc 3 (?7i - 773)) . 

m 

3.6 Second order curvature perturbation and summary 

So far, we have only discussed the linear part of the curvature perturbation since it is the only 
term that picks up contributions from cubic vertices at tree-level. The second and higher 
order terms in £ also contribute to the bispectrum, however, through the combinations such 
as (C(i)C(i)C(2)) ■ Since it is impossible to examine at all orders if they give any contribution 
within the order X 2 , here we just look at the second order term and check that they do not 
become dominant over the 3-vertex contributions derived in the previous subsection. 

First, we note that ignoring higher order corrections in eH,ip,VH,ip and the terms with 
spatial derivatives and using (2.17), W6 Cciii rewrite as 

4 C(d - |Mi)C(d) <% - 1 (C(i),*s« + CmjbV) . 
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Comparing eq.(3.3) with 



i 61 3 o f q 2 
7T H 7T i=Hrf a 



we see the second term is suppressed by a factor of e#. Throwing it away, eq.(2.18) yields 

-C(2) - ^<V>(2) + 2— ^C(i) - 2C(i)- 



Expanding the energy-momentum tensor up to second order, we find 
2 

P + P 



Sp {2) - -2T ( ° 2)0 + — -T^T^ 



V , 2c 2 (/% 1 6c (/ 2 ),„ / , 12c / , / 2 , , 2 > 

The appearance of (f)(2) forces us to look into the constraint equations at the second order. 
In fact, they are not too bad for the scalar perturbations in the flat gauge. The relevant 
equation is obtained from variation of Ni in the ADM formalism and reads 

ifj = -vVj + <P,iB,i<p d + f 2 L?F a ). 

Expanding it to the second order, we find 

2W0( 2 )j = (Wd) + V 2 ^ 1 )) - (3-20) 
f 2 

+ A; + ^2 ( 2a ' a J + T >"« + r J v2r ) • 

In the end, its contribution is subdominant. Keeping the leading order terms in slow-roll and 
discarding higher spatial derivatives, we obtain 

C(2) = ^vr' 2 + J- (24ZV - YlHX^ol + 3#V«' 2 ) + 2C ( 2 1) - 27^(1)- (3-21) 
The contribution to the bispectrum is 

(C(x)C(y)C(z)) ~ (C(l)(x)C(l)(y)C(2)(z)) + (2 perms). 

The first and the last terms in (3.21) are subdominant. The term quadratic in £n) gives 
contributions such as 

(C(i)(x)C(i)(z))(C ( i)(y)C(i)(z)) ; 

which exist regardless of the dynamics and give |/jvx| i$ 1- The rest are the generic effects of 
the background gauge fields. Looking at (2.19), we see that the leading order contributions 
in X is quadratic, which involve terms such as 



f 3 



X 2 2( 7 r(x) 7 r(z))(7r(y)7r(z)) + -#VV(x)a'(z)) (a'(y)a'(z)) 
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and 

-JV^W ((7r(x)7r(z))(a'(y)«'(z)) + (^(y)vr(z))(a'(x)a'(z))) . 

At the leading order, tt and a' are essentially just Uk{rj) and /^(ry)/?? 4 , therefore their contri- 
bution will be constant of |/jvi| ~ 0(I 2 ). 

In the end, our result is summarised as follows. At the order of Z 2 , the tree-level 
amplitude of the three point function in the super horizon limit becomes 

3e H 4 I 2 

{0fci0fc 2 0fc 3 ) 9 ^TeTeTe ( fc 3^3a + A^fe + ^l^c + 2 perms) 

+ + $ (ln(-^)) 3 

where if is a reference momentum, say K = \{k\ + &2 + ^3). While the ambiguity of K 
arising from the lower limits of the integrations leads to errors of order ln(ki/kj),i 7^ j, 
for the wavelengths of interests, this should be of order 10. Since there are many other 
contributions of similar order which we have already ignored, it does not make sense to 
overly worry about this reference momentum. As we can see, the bispectrum is of local 
shape. In order to estimate the f^L in the squeezed limit, we quote the result from ref.[61] 
for the power-spectrum 

{CD "f^ (l + 18V6X 2 (H-k V )) 2 ) . (3.22) 

Under the condition 

-I\n{-Krj) <C 1, (3.23) 
which implies we can replace e v with en, we obtain 

f NL ~ 972X 2 jV|, N K = - In (-K V ) . (3.24) 

This result qualitatively agrees with the one derived in ref.[52]. One cannot drop the re- 
striction (3.23) since we are sticking to the leading order contribution in X. Noting Nk 
is the e- folding number experienced by the relevant modes after horizon exit, it has to be 
greater than 50. Tentatively setting Nk ~ 100, the current observational bound Jnl^S 100 
translates to 

I 2 < HP 7 . 

This is a strong constraint. The physical interpretation is that the energy density of the 
background gauge fields is seven orders of magnitude smaller than the kinetic energy of the 
inflaton, which is already suppressed by e v with respect to the vacuum energy. 

However, there are a few unsatisfactory features in this result. The first is the limitation 
(3.23) that itself already restricts our discussion within X 2 <C 10~ 4 . Furthermore, the fact 
that may grow indefinitely as long as inflation continues sounds unpleasant considering 
the classical stability of the quasi-de-Sitter background. It is distinct from the infrared 
divergence discussed in ref. [52] which concerns the back reaction of the quantum fluctuations 
and loop corrections which is beyond the scope of the present article. The divergence is 
already there at the tree-level calculation. Motivated by this, in the next section, we shall 
give a more careful analysis on the superhorizon dynamics of the fluctuations. 



- 22 - 



4 Non-perturbative treatment of the quadratic vertices 



It is clear that the above approach based on the perturbative expansion in terms of X bares 
a limited applicability even if X <C 1. From the point of view of the classical stability of 
this inflationary regime shown in ref. [53], the apparent indefinite growth of the correlation 
functions after horizon exit should halt sooner or later if all the relevant effects are taken 
into account. In the Lagrangian (3.1), we have regarded the quadratic interaction terms 
as perturbative corrections along with the cubic ones. In this way, the proper tree-level 
amplitude involves an infinite number of Feynman diagrams generated by those quadratic 
vertices. While we have avoided this issue by focusing on the leading order contribution in X, 
one should expect a convergent result if the higher-order corrections are treated appropriately. 
For this purpose, we investigate the linear perturbation more closely and show that both the 
power spectrum and the bispectrum become constant in the limit of r] — > despite the 
appearance of logarithmic divergence ln(— krj) in the perturbative analysis. 

4.1 Linear evolution equations and their superhorizon solutions 

The equations of motion at linear order are given by 

1 , „ l 2 2 , 2 , 12X 2 6y/2X , 
Jp^2 K + k n k) ~ Jp^^k ~ Jp^^ = ^-a fc (4.1) 

(3r ? V fc )' + 2fcV« fc = ^<- (4.2) 

It turns out that one can write down analytic expressions for the solutions in the superhorizon 
limit. Ignoring the spatial gradients, the second immediately integrates to give 



, c , 2v/2X 



where cq is an integration constant. We suppress its ^-dependence since there should be 
no confusion as far as the linear theory is concerned. The same applies to the rest of the 
integration constants. Plugging this into the first equation, we derive 

„ 2 , 12X 2 GVZHIcq 

7Tj. TT k -\ 5— 7Tfc = k 

whose general solution can be written as 

nk = "Tlx (co + c+ { ~ kv)P+ + c " (_A;r?r) (43) 

with two arbitrary constant c±. The power exponents are given by 

3 ± V9 - 48X 2 
P± = 2 • 

The corresponding a is 



+ 2(-Af (d + ^ (-k V y p - + — (-k V y p+ ) (4.4) 

V P- p+ J 

with the fourth integration constant c\. We used the relations 

p± - 3 = -p T . 
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4.2 Canonical mode functions 

When the off-diagonal terms in the quadratic Lagrangian are taken into account, the intro- 
duction of mode functions is not so straightforward as with independent free fields. In this 
subsection, we look into the canonical formulation of the field theory. From the Lagrangian 
(3.1), we read off 

7T = air and a = V3fa 
as the canonically normalised field variables. Their conjugate momenta are given by 

Pit = 7r', 

p a = a -\ TT 

V 

and we impose the canonical commutation relations 

[7r(T,x),p 7r (T,y)] =i6(x.-y), [d(r, x), p a (r, y)] = i5(x-y) 

with all the cross commutators being zero. To diagonalize the Hamiltonian, we introduce the 
creation and annihilation operators 



^a^(p-q), 0,6 = 1,2 
and expand the field operators in terms of the mode functions: 

d 3 k 



7t(? ? ,x)= J (ugC^a^e** + 4*^)0^6 



ikx 



a=l,2 ' 



Here, (ui,vi),a = 1,2 are two independent solutions of the eqs.(4.1) and (4.2). As an 
example, for the superhorizon solutions derived in the previous subsection, the mode functions 
become 

(eg + c\ (—kr]) p+ + ct (—krj) p ~) , (4.5) 



7/ V 6 P- P+ J 

which are characterised by eight complex constants. In this way, we see that each field 
operator may excite two different particles. Conversely, for each particle species a, there 
are two associated mode functions u^, and vi. This simply reflects the fact that the fields 
themselves do not define particles when there is a quadratic mixing term. This formulation 
is consistent as long as the mode functions satisfy the following conditions arising from 
the canonical commutators (from here on, the summation convention for indices a,b, - ■ ■ is 
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assumed) : 

u a u a *' - u a *u al = i, 

v a v a*> _ v a* v a, + ra^* _ = ■ 

rj 

u a v a * - u a *v a = 0, 

U <V*' _ U ** V *> = 0, 

y a u a*> _ v a* u a, = ^ 

u a> v a*, _ u a*, v a, + < ^f$Z , al ^ _ = Q _ 

It can be checked that they are preserved by the evolution equations (4.1) and (4.2) if they 
are satisfied at an initial time. 

Later it proves to be useful to write down these conditions specifically for the superhori- 
zon mode functions. Those six equations translate into algebraic conditions on the integration 
constants: 



2X 2 



(4.7) 



+ + A: 3 (p 

44* - eg* 4 = c?4* - cf 4 = 0. (4.9) 
4.3 Matching with the de-Sitter mode functions 

From the solutions (4.3) and (4.4), it is almost obvious that the power spectrum should 
converge to a constant proportional to |cg| 2 - In order to estimate its magnitude, however, 
one needs to determine the constants 4 , a = 0, 1, ± from appropriately initial conditions set 
deep inside the horizon. As a first approximation, we can match the de-Sitter mode functions 



\/2k V k V 



^2 ( 1 l _ \ e -ic 



y/2c s k V c s kr] 

with the superhorizon counterparts (4.5) and (4.6) at the horizon crossing krj = —1. We 
evaluate them and their first derivatives and equate each other. This leads to the following 
eight equations: 

Ti) a 

c a Q + c a + + c a _ = --^(l + i)e i , 

c a a 4 c a _ 5 a 2 / i 

6 p_ p+ ^/2Ac s k 6 V °s 

-4 + (p + - 1) 4 + (p_ - 1) c a _ = -j±e\ 



-| + 24 + g±^4 + g— ±<? = /° 2 ,(l + ^i)e 
6 1 p- + P+ ~ V24cTF V 2 
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These easily solve as 
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4 



4 



21 



(1 + 



1 



(6 + 7c s i) 



k 3 Pt ~ P± 



3Vk 3 



l + < + 



37 



1 + i)p 



4 



p+p 

(3 + 4c s i) 



P± ~Pt 



4 



We can now estimate the final amplitude of the two point function. Clearly, the domi- 
nant contribution at late time comes from Cq's. The mode functions asymptotically approach 



v a k (v) 



Using 



we derive 



Hr] 2 
3e H 



4 + 6X 2 



V 



-7T + 
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—a 

V 



44^ 



\co\ 



'1+^ 2 + 



103 



144c,X 2 



(4.10) 



Although the exact numerical factors should not be trusted due to the errors coming from 
the matching, the dependence on I 2 is generic (we will confirm this numerically). The fact 
that the final amplitude is divergent in the limit small I 2 might look worrying. But it is 
also the case in the standard single-scalar inflation where the power spectrum is formally 
infinite when the background is exactly de-Sitter. The same nature can be seen taking the 
limit of 1 — > and large e-folding number (— kr] — > 0), rewriting the amplitude of the power 
spectrum as 

(CD - O(l)^ ~ 7—^- (4-11) 



In the last approximation, we used (2.5) and 



X< 1 ^ — ~ 1. 

This result is reasonable: the dependence of the power spectrum on the parameters is the 
same as the single-scalar inflation except that is now replaced by ejj — Recalling that 
it represents the energy density of the gauge fields in the unit of H 2 , the power spectrum 
is inversely proportional to the energy density of the background gauge fields instead of 
the background scalar kinetic energy. We emphasise that the only approximation we used 
to derive (4.10) was the matching with de-Sitter mode functions. Hence, we expect the 
expression is valid even for 1 > 1 in the super horizon limit. 

Comparing the expression (4.10) with the perturbative result (3.22), one can estimate 
the time when the power spectrum settles down to a constant value after horizon exit. We 
simply equate these two in the limit of small 1 and infer that 



(4.12) 



Beyond this point, (g!) is conserved as it is in the usual adiabatic perturbation. 
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4.4 A proof of convergence of the bispectrum 

Under the general conditions (4.7) - (4.9) arising from the requirement of canonical commu- 
tation relations, one can show that the tree- level amplitude of the three-point function is 
convergent in the superhorizon limit too. 

First of all, we introduce the mode functions for £ by 

d 3 k 



Cfa,x) = J |A( a(7?)e *-x a a +cr(r?)e -ik.x a at 



where 



C k = [3 (i + 1 2 ) eg + {p + + 3) 4 (-k v f+ + (p_ + 3) e (-k^y-] . 

Since there is no quadratic mixing term (that are already taken into account in defining the 
mode functions), by taking £ itself as an independent variable, the relevant vertex must take 
the form 



f C 3 
H 1 } = const / d 3 w-± 



Then, the tree-level amplitude becomes 

(Cf) = eonst J" ^9 (c?Wft(»7)Cg(i7)Ct(»Zi)4(%)C?,(m)) • 



Notice that 
3fif 2 T 2 

^||^cr(«(m) = 9(i+x 2 )|csi 2 

+ {p + + 3) 2 \c a + \ 2 (-kv) p+ (-kvi) p+ + (P- +3) 2 \c a _\ 2 (-kv) p -(-kvi) p ~ 

+3(1 +i 2 )( P+ + 3) [44*(-M p+ + cg*4(-W + ] 
+3(1 +i 2 )( P - + 3) [c$cT(-k V1 y- + ^cti-kvy-) 

+ ( P+ + 3)(p_ + 3) [c a + c^(-kr,)P+(-k m y- + 4*c^(-fcr ? ) p -(-fcr ?1 ) p +] . 
We can see that the lowest power of the integrand must come from 

^ {Q:(V)QM)) » (ctMd 2 (m)) K {CZlWClAm)) + 2 perms. 
The time dependence of its dominant contribution is given by 

m A ^(Cr(vKk(m)) oc 3(1 +x 2 )( P+ + 3) (44* - 4*4) - (-A^rhr 4 

+3(1 +X 2 )(p„ + 3) (4c»* - 4*c^) [(-fcr/xf- - (-^i) p i r?r 4 

+(P + + 3)(p_ + 3) (c a + c a * - c a + *c a _) l(-k V y+ (-k Vl y- - (-krif- (-k m y+] r/7 4 . 

However, the first and the second lines vanish because of the conditions (4.9). The total 
power of the remaining term is —1 and it contains a positive power of rj, which implies the 
integration in the limit —rj — > is convergent. 

Now the question is to evaluate the magnitude of this asymptotic constant value of 
the three point function. Although not entirely trustworthy, it is reasonable to expect that 
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the bispectrum becomes constant around the same time as the power spectrum does. As a 
rough estimate, we substitute (4.12) into the perturbative amplitude (3.24), and find that 
the final value of Jnl is expected to be independent of I. The numerical factor can easily 
vary between 0(1) and 0(100) depending on the way we make the approximations. In the 
next section, we use a numerical method to derive its dependence on X and constrain it in 
the regime where the perturbative result (3.24) cannot be trusted. 

5 Numerical calculation of exact tree-level amplitude 

Following from the previous section, here we treat the quadratic vertices non-pertubatively 
with the only difference being that now we calculate most of the contributions numerically. 
The aim is to negate the need for making any approximations and therefore to make our result 
more quantatively accurate. In the analytic results from section 4, we derived the qualitative 
features of power spectrum and bispectrum by assuming that at horizon crossing the mode 
functions are those of the free de-Sitter case and applying the superhorizon approximation 
(krj = 0) as soon as the mode crosses the horizon (krj = — 1). We have been able to estimate 
the final amplitude of power spectrum, the time of transition from the perturbative regime 
discussed in section 3 to the one dictated by the superhorizon mode functions, and prove the 
convergence of the bispectrum in the limit krj —> 0. Nevertheless, neither have we obtained 
a reliable estimate of the convergent value of bispectrum nor the behaviour between horizon 
crossing and convergence; this is clearly unsatisfactory in the light of precision cosmology. 
Here we calculate the exact mode function, first setting the tt and a fields in the Bunch-Davies 
vacuum deep inside the horizon, solve the coupled linear equations of motion numerically until 
the modes are far into the superhorizon regime. At this point we switch to using the super 
horizon equations of motion and use the analytic solution - this is simply to avoid numerical 
instabilities encountered in this calculation. Now we only use the analytic super horizon 
solution for — krj <C 1, so the error introduced by doing so is negligible. 

Ultimately, we will be interested in the value of /jvx (squeezed) here. Factors of JqH 
present in the Lagrangian (3.1) will be absorbed into the definition of the a field, and the 
overall multiplicative factor of H~ 2 in front of the Lagrangian will not affect the value of 
Inl- We therefore set H = 1; for quantities such as the power spectrum or bispectrum, 
reintroducing will H will be a matter of an overall multiplicative factor which will be included 
in the plots. When reintroducing H, I will need to be replaced with X/H. 

This leaves the factors of en and in the Lagrangian and the definition of the curvature 

perturbation; in the numerical calculation below, they will be set to 1. It can be easily 

seen that these two parameters can be reintroduced at the end as an overall multiplicative 

factor of en/^ifi = 1 + X 2 for the value of /nl computed. The £ modefunctions, power 

_i 

spectrum and bispectrum will need to be multiplied by H(l + I 2 )~ 1 e v 2 , H 2 (l + X 2 )~ 2 e~ 1 
and H 4 (l + X 2 ) _3 e~ 2 respectively to restore the dependence on these constants. 

5.1 Sub horizon linear evolution and initial conditions 

Again, the linear equations of motion (4.1) and (4.2) are used, this time keeping the gradient 
terms. Since the evolution equations do not admit an analytic solution, we will find solutions 
numerically. For computational convenience, the canonical variables used in this section are 
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Figure 1. Plots of the modcfunctions during subhorizon, with kr\ along the x-axis. For all plots 
on this page, the solid line represents the real part and the dashed line represents the imaginary 
part. The different plots are for different values of I; from left to right: 2" = 0.1, 0.5, 1 and 10 
respectively. For smaller values of I, one can observe the characteristic oscillation and its damping 
towards horizon exit of de-Sitter mode functions, while the behaviour near the horizon is significantly 



different for X = 10. 




Figure 2. Plots of the Q, modefunctions during subhorizon, with kr\ along the x-axis. The subhorizon 
dynamics appears to be similar between Q. and Q.. 



Figure 3. Plots of the Q, modcfunctions during super horizon, with cosmic time (— ln(— fa])) along 
the £-axis. While the evolution of individual mode functions significantly depends on the value of X, 
all of them settle down to constant in agreement with the analytical results. 



Figure 4. Plots of the (| modefunctions during subhorizon, with cosmic time along the x-axis. On 
superhorizon scales, Q and £| evolve differently. 

7r and a, and so their conjugate momenta are given by 

Pw = aV , 

p a = 3a~ A a' — 6V21ir . 
The initial conditions for the modefunctions are given by the Bunch-Davies condition, 
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expressed here in terms of the canonical variables for each k mode as: 



= —!==(» — kr})e~ ikv , p% = %a? J %i V e- ikr > , (5.1) 



^(cskrj^-i^e-^, (5.2) 



p% = ai J— (-ic 2 s k 2 r,- 2 - 3c s k V - 3 - 3ir)- 4 )e- ic ° k « . (5.3) 

The point here is that these are the conditions required on the mode functions for the fields to 
be in the Bunch-Davies vacuum deep inside the horizon, and for the canonical commutation 
relations to hold. Given the definitions of the conjugate momenta, the initial conditions for 
solving the linear evolution equations will then be given by (5.1) and (5.2) along with 

a% = 2 3 3 (-ic 2 s k 2 V - 2 - 3c s kr]- 3 - ^'^e'^ + ^7=("^~ 3 + ir]- 4 )e~ ik ^ (5.4) 

in place of (5.3) for some — krj S> 1. For the results in this section, the initial conditions for 
the modefunctions were set at (r))imt = ~~ 1000. 



Figure 5. Time evolution for the power spectrum (Cfc 2 )- The first line shows the time evolution 
until the power spectrum becomes constant and the second line shows the time evolution for a some 
time after horizon crossing, both plotted against cosmic time. The second line provides a comparison 
between the perturbative result from ref. [53] and the numerical results from section 5 for the time 
evolution of (Ck 2 ), with the numerical results in solid lines and the analytic ones in dashed. From left 
to right: 1 = 0.01, 0.1, 0.5 and 1. While the analytical results fit the initial growth fairly well for 
X = 0.01,0.1, the discrepancy is significant for larger I or late time. 



5.2 Numerical calculation of the £ power spectrum 

Now we calculate the power spectrum for the curvature perturbation. A similar analysis has 
been carried out in ref. [49], so the results in this section are to recap these results, and to 
verify that these results are consistent with the perturbative expression for the curvature 
power spectrum[53]. Using the tt and a mode functions we are able to define the £ mode 
function as 

/-a V ( 'aJ , 3(1 + 21 2 ) 3 3„aA /c c\ 

Ck = "sTi V k + — v — " fc " 7^ 1 k ) ' ( 5 
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(t-Log[kl) I 2 



Figure 6. Time evolution of (Ck 2 ) rescaled according to analytic expectation; the (overlapping) 
dotted and dashed lines are X = 0.001 and 0.01, and the solid line is I = 0.1. The time coordinate 
has been rescaled by I 2 and the power spectrum amplitude by I 2 . It clearly indicates the time of 
transition to the constant regime — ln(— kr)) ~ I~ 2 . 

final P f 

1(1" 
10 s 
10" 
1000 
100 



0.005 0.010 0.050 0.100 0.500 1.000 

Figure 7. Values for the power spectrum ((|) after settling down to a constant, plotted as a function 
of X, The solid line is the numerical result, and the dotted line is the analytic expression from section 
4. There is a good agreement for I < 0.2. 




and hence 

(Ck 2 )=C k C k . (5.6) 

The C modefunctions are plotted in Figures 1, 2, 3 and 4, while the numerical and analytic 
results for the power spectrum are shown in Figures 5, 6 and 7. The perturbative solution 
for the time evolution of (Cfc 2 ) is shown to be useful only for small X (Figure 5), while the 
analytic estimate for the final value of the power spectrum, derived in the section 4, is valid 
for values of X up to around 0.2 (Figure 7). The lack of quantitative agreement beyond 
X ~ 0.2 is presumably due to the error arising from the matching since for larger values of X, 
the numerical calculations (Figrues 1, 2, 3, 4) show a significant deviation from the de-Sitter 
mode functions around horizon crossing. The characteristic time scale for the time evolution 
for (Cfc 2 ) before it reaches constant is shown to be X -2 for X < 0.1, in agreement with the 
analytical estimate (4.12) from the previous section (Figure 6). 

5.3 Numerical calculation of the ( bispectrum 

By solving the coupled linear evolution equations we in effect include the contribution from 
the infinitely many tree-level Feynman diagrams coming from the quadratic H q term, and 
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hence obtain a result correct to all orders in X (provided loop contributions are negligible). 
Therefore, the exact tree-level amplitude for the bispectrum, by standard application of 
Wick's theorem, is given by 

(c kl cM = - i" ^ (c fe o r(r?)ci:wa:(^x(^X(^)<3(m)) 

+ 48V2X fdm* (arWCSWCg^X^iX^iX^i) + 2 perms) 
- 12x/2 |"d7 ?1 r ? ^(c^( ?? )C^(^)a;(^)4 1 (m)<(m)<(m) + 2perms) . 

In particular, we now only have to compute 1-vertex terms. 

The evaluation of the integrand turns out to be a numerically unstable process suffi- 
ciently far outside the horizon, requiring a very precise cancellation of terms. It therefore 
becomes impractical to carry out the calculation with the numerically solved mode functions 
beyond a certain point. To overcome this difficulty, for —krj < 10~ 5 we switch to using the 
analytic superhorizon solution discussed in the previous section. The only difference is in 
the matching of the analytic superhorizon solution; here we evaluate the numerical mode 
functions (and their time derivative) at — krj = 10~ 5 and use these as the matching condi- 
tions for the analytic superhorizon solution. Then, for —krj < 10~ 5 the time integrals in 
the above expression are computed analytically, therefore avoiding the problem of numerical 
instabilities. When performing the first stage of this computation (the numerical stage), we 
employ the technique recently developed in [64]. As we will see later, the bispectrum (or more 
precisely, the shape function) is peaked in the squeezed limit, and therefore to concentrate 
on the salient features we will restrict most of our analysis to the squeezed limit. 

We start by crosschecking our numerical calculations against the perturbative results 
from section 3 (Figure 8), for small X (= 10 -3 ) for some time after horizon exit where the 
perturbative treatment of the quadratic vertex H q is justified. This underwrites the overall 
consistency between the analytical and numerical methods. 

In Figure 9, we confirm the convergence of the bispectrum generated by each cubic 
vertex. As one can see, while H c is the dominant contribution in the perturbative regime 
as it grows fastest (ln(— krj)) 3 , it is overtaken by H B when the perturbative approximation 
breaks down. It also exhibits an approximate X -4 scaling of the final value of bispectrum, 
which is equivalent to the X independence of f^L that was inferred at the end of the previous 
section. The characteristic time scale for the transition to constant is again shown to be X -2 . 

From the (ln(— krj)) 3 perturbative growth in the bispectrum, one may expect that the 
superhorizon contribution to the bispectrum dominates over the subhorizon contribution; in 
Figure 10 we verify that this is indeed the case for most values of X. The bispectrum is 
evaluated at r/ = 0, and the krj < — 1 (subhorizon) and krj > — 1 (superhorizon) contributions 
to the time integral are plotted separately as a function of X. The two become comparable 
only as X reaches order unity, and for X < 0.1 the subhorizon contribution is negligible. 

In Figure 11, we plot the intermediate time evolution of Jnl, with the numerical cal- 
culation on the left panel and perturbative result on the right. For <C X -1 , the power 
spectrum is essentially constant and /jvx grows as N%. When 7V& > 0.1X -1 , the power spec- 
trum starts to be overtaken by the correction term and scale as N?, which results in the 
peak around iV& ~ 0.3X -1 . The maximum value appears to scale as X" 1 , which means it 
may well be observable for a small value of X. Since these peaks occur on time scales oc X -1 , 
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Figure 8. Comparison between the analytic (perturbative) and numerical values of the bispectrum 
in the squeezed limit, for X = 0.001 just after horizon crossing, showing very good agreement; from 
left to right are the contributions from H A 1 H B and H c respectively, with the numerical value in 
solid lines and analytic in dashed. For all the results in this section, the squeezed limit is evaulatcd 
by taking the bispectrum in the configuation kt x 10 2 = k% = ^*3- 



Figure 9. Late time evolution of the bispectrum in the squeezed limit. There are two properties here 
which cannot be seen from the perturbative calculation; the contribution to the bispectrum from the 
second vertex becomes negative, and the contribution from all three vertices become constant after 
some number of e- folds ~ Z~ 2 . For small I, the final value of the squeezed bispectrum exhibits a I~ A 



the time dependence (and therefore scale dependence) of Jml around this maximum can be 
understood by the perturbative results where analytical expressions are available. 

However, for later times > X~ 2 this is no longer the case, with /tvl turning negative; this 
behaviour is shown in Figure 12, where we plot the late time evolution of Jnl m the squeezed 
limit for different values of X. The final convergent value appears to be —5/3 3 , independent 
of T. It is also seen in Figure 13 where the final value is presented as a function of X. We 
suspect the cause of irregular behaviour for X > 0.1 is due to the significant contribution 
from the subhorizon evolution. We also note that given that there are two scalar degrees of 
freedom here, the single field consistency relation [2] does not hold. 

We conclude this section by mentioning a few words about the shape of the bispectrum. 
As can be seen from Figure 13, the bispectrum is peaked in the squeezed limit, which is 
expected given the fact that it is predominantly determined by the superhorizon evolution 
which tends to generate local bispectra. Provided we wait until all relevant modes have 
become constant (as is done for the plot), it is perfectly scale invariant too. This can be 
understood by noting that since the background geometry is de-Sitter and the interaction 
terms are de-Sitter invariant, the correleation functions for the perturbations which have 
become constant are scale invariant; in particular, for the modes which have settled down 
to the final value, both the power spectrum and bispectrum are scale invariant. By similar 
arguments, where we have plotted the time evolution of any quantity such as Jnl, they can 

3 This peculiar fractional value is due to the | term in the definition of /jvl, and without this the final 
value would have been -2. 




scaling. 




Figure 10. Relative importance of the subhorizon and superhorizon contributions to the final value of 
the squeezed bispectrum. As X approaches 1 , the subhorizon and superhorizon contributions become 
comparable. For the plots of the H B interaction and combined interactions, the negative of the 
bispectrum was taken for the purposes of taking a log plot. We were unable to obtain reliable results 
for the subhorizon contribution for I > 0.1, so they were not included in the above plots. 




Figure 11. Intermediate time evolution of Jnl in the squeezed limit for I = 0.001,0.003,0.01 and 
0.1. Jnl first grows as 7V| and eventually peaks around Nk ~ 0.3/1. The peak value appears to 
be ~ I -1 . Then it monotonically decays until eventually settling down to a negative constant. The 
behaviour of the peak, at least for I < 0.01, can be understood perturbatively since the timescalc for 
the peaking of /nl is smaller than that of the breakdown of perturbation theory. 



be used to read off the scale dependence at any given time. 
6 Concluding remarks 

We have studied the perturbation of a model of inflation where a stable isotropic phase of 
inflation is realised by a scalar field coupled with a triplet of Abelian gauge fields. We derived 
the general action for scalar perturbation up to cubic order and identified all the relevant 
terms in the limit of vanishing slow-roll parameters. Using the standard method of in-in for- 
malism, we first treated both the quadratic and cubic vertices perturbatively and computed 
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Figure 12. Late time evolution of Jml in the squeezed limit for I = 0.001,0.06,0.1 and 0.2. The 
x-axis is e- folding number after horizon crossing rescaled by I 2 , showing the characteristic timescale 
for fjsiL to become constant. The final value is independent of X. 



Final squeezed f M 





Figure 13. Left: the final value of Jnl in the squeezed limit, at least for 1 < 0.1, is independent of 
1 appearing to take the value — |. Right: The shape of bispectrum for 1 = 0.01, final value. What 
is plotted here are contours for (fcirC2rC3) 2 i?(fci, fc 2 , fc 3 ), with the three axes being fci, and k^. It is 
peaked in the squeezed limit and is scale invariant. 



the bispectrum at the leading order in the expansion parameter I. The resulting expression 
was consistent with the previous studies and f^L in the squeezed limit was shown to be 
proportional to 1. 2 N^. where JV^ is the e-folding number after the relevant modes exit the 
horizon. We then pointed out the limited applicability of this approach even for X 1 and 
rectified it by introducing the exact linear mode functions which take into account the effect 
of the infinite number of tree-level diagrams generated by the quadratic vertices. Solving the 
linear evolution equations analytically in the super horizon limit, we proved that both the 
power spectrum and bispectrum are convergent in the limit Nk — > oo and estimated their 
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dependence on X by the means of matching with the de-Sitter mode functions at horizon 
crossing. In order to obtain a more quantitative estimate of the bispectrum and /nl, we car- 
ried out an extensive numerical analysis employing in part the recently developed technique 
[64]. We confirmed the analytical results and found a number of interesting features. The 
peak value of /nl in the squeezed limit appeared to scale as X~ l and may well be observable 
in the forthcoming observations. The convergent value of f^L m the limit of Nk — > oo was 
shown to be negative and independent of I. The shape of bispectrum is predominantly local 
as is expected by the fact that it is dominated by the superhorizon contribution. 

We first emphasise that our analysis here is complete at tree-level. It indicates the 
overall consistency of this model in the classical regime. Regardless of the value of I, the 
model can be consistent with the current observational constraint since the final value of /jyx 
is of order unity while it offers a realistic possibility of detecting a significant non-Gaussianity 
if we happen to catch the period shortly after the horizon crossing in our observable window. 
We note that the smaller the value of I, the greater the peak value of /n l and the longer it 
takes for /jvl to settle down to the constant. 

One can expect that the same qualitative features will also apply to the anisotropic 
models where the background is permeated by one or two gauge fields with non-vanishing 
vacuum expectation values. The difference is that there the vector fields also contribute to 
the spatial anisotropy and there is a strict upper bound for I. It is going to be difficult to 
repeat our analysis for anisotropic models since the consistent perturbative expansion requires 
inclusion of vector and tensor modes which are coupled to the scalars through the background 
anisotropy. For this purpose, it would be interesting to look into the relation between our 
results and the delta- N formalism [65]. In fact, the isotropic case can be regarded as a 
particular two-scalar model and the formalism should apply without any problem. Since the 
convergence of the power spectrum and bispectrum is based on the superhorizon evolution, 
the delta- N formalism will reproduce them in a more elegant manner. Since its mathematical 
basis resides in the equivalence of the superhorizon curvature perturbation to the background 
evolution of the FLRW universe [66], an appropriate extension to anisotropic backgrounds 
sounds plausible and can be a powerful tool to handle the complicated interactions among 
different modes. 

Another important theoretical issue is consistency of the quantum field theory in the 
existence of background gauge fields. The authors of [52] claimed that the infrared contribu- 
tion of the one-loop diagrams can be interpreted as the rescaling of the background vacuum 
expectation value of the gauge fields so as to take into account the quantum mechanically 
created modes that froze in outside the horizon. Although we have not discussed this issue 
in the present article, it will be certainly an interesting direction of further research. 

Finally, it is in principle straightforward to extend our analysis to inflationary mod- 
els with non-Abelian gauge fields, either the one based on gauge-kinetic coupling [44, 54] 
or Chern-Simons coupling [57]. Given the qualitative similarity between Abelian and non- 
Abelian models when all the vertices are treated perturbatively, it is natural to expect that 
a similar convergent result in the limit of large e-folding can be established, although it is 
mathematically far from obvious. From a phenomenological point of view, it would be impor- 
tant to clarify the difference among different scenarios so that one is able to observationally 
distinguish between them. 
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